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ABSTRACT 

o, 

^ , We present a galaxy evolution model which incorporates a physically mo- 

^ ! tivated implementation of AGN feedback. Intermittent jets inflate cocoons of 

^! radio plasma which then expand supersonically, shock heating the ambient gas. 

2 '. The model reproduces observed star formation histories to the highest redshifts 

c/3 , for which reliable data exists, as well as the observed galaxy colour bimodality. 

Intermittent radio source feedback also naturally provides a way of keeping the 
black hole and spheroid growth in step. We find possible evidence for a top-heavy 
^ \ Initial Mass Function (IMF) for 2; > 2, consistent with observations of element 

\ abundances, and sub-mm and Lyman break galaxy counts. 
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Subject headings: galaxies: evolution — galaxies: formation — galaxies: active 



1. Introduction 



X 

I A remarkably simple and self-consistent picture of structure formation via growth of 

the gravitational instability has emerged over the past few decades. Recent measurements 
of the cosmic microwave background (CMB) power spectrum by the Wilkinson Microwave 
Anisotropy Probe (WMAP; Spergel et al. 2007), and the galaxy clustering power spectrum 
observations from the two degree Field Galaxy Redshift Survey (2dFGRS; Colless et al. 
2001) and the Sloan Digital Sky Survey (SDSS; York et al. 2000) have confirmed the cold 
dark matter model with a non-zero cosmological constant (the ACDM model) as providing 
the current state-of-the art description of structure formation. This cosmogony is consistent 
with the rate of expansion of the Universe as inferred from supernovae observations, present 
mass budget of the Universe, observed baryonic fraction in rich clusters, structure of the 
Lyman-a forest at z = 3, and the model of Big Bang nucleosynthesis (Spergel et al. 2007). 
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Structures grow via amplification of over densities. Gunn & Gott (1972) were the first 
to use the spherical top-hat model to track the growth of clusters. Press & Schechter (1974) 
extended this analysis to growth of Gaussian fiuctuations by smoothing the field on dif- 
ferent scales to predict the dark matter halo mass function. The results turned out to be 
in remarkable agreement with later numerical simulations that followed the evolution of 
structure subject to the gravitational instability (e.g. Lacey & Cole 1994; Somerville et al. 
2000). Sheth et al. (2001) refined the Press-Schechter formalism (to the so-called Extended 
Press-Schechter formalism, EPS) to consider ellipsoidal, rather than spherical, collapse. The 
merger histories of the dark matter haloes can be extracted by employing either percolation 
algorithms or Monte-Carlo simulations. With the merger trees in place, it makes sense to 
talk about individual halo accretion histories. Analytical fits to halo mass functions (e.g. 
Jenkins et al. 2001) and average mass accretion histories of individual haloes of a given 
mass at the present epoch (e.g. van den Bosch 2002) are now available. In addition, the 
density and temperature structure of dark matter haloes has been studied extensively using 
numerical simulations. 

Observable galaxy properties are due to the baryonic component of these structures. 
White & Rees (1978) suggested that haloes accrete baryonic matter along with the dark 
matter. The baryons are shock heated to the virial temperature of the halo, T^ir, and 
thermalization of dark matter particles supports the halo from further collapse. As the hot 
gas cools, it is no longer supported by radiative pressure, and sinks to the centre of the halo 
potential. Once cold, the gas suffers from the gravitational instability, and collapses to form 
stars (Kennicutt 1989). 

While very attractive, this picture is at odds with the observable galaxy properties. 
In the hierarchical assembly paradigm lower mass galaxies form first, and it is the most 
massive galaxies that should be undergoing vigorous star formation at the present epoch. 
However, observations (e.g. Silk 2002 and references therein) show these massive galaxies to 
be systematically redder (as calculated from their optical colours), implying stars in these 
systems are older, than their less massive coun terparts. Infra-red colo urs of these objects are 



consistent with passive evolution since 2; ~ 1 (IBender fc Saglialll999l ). Using the rest-frame 



UV fiux as a p roxy for instantaneous star formation rate, derived star formation rates peak 



downsizing" scenario (jCowie et al. 



at z ~ 1 — 2 (IGiavalisco et al.l 120041). c onsistent with this picture. Thus, in this "cosmic 



19961 ) the higher mass galaxies are assembled first. 



Further evidence for suppression of star formation in massive galaxies comes from ob- 
servations of the local optical luminosity function. Under the assumption of a constant 
mass-to-light ratio, the number counts of the most massive and least massive galaxies are 
significantly overpredicted when compared with observations (Baugh 2006; Eke et al. 2006). 
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While reionization and supernovae feedback provide sufficient heating to reconcile counts at 
the faint end of the luminosity function, this is insufficient in the most massive galaxies. 

The overprediction of stellar content in the most massive structures at the present epoch 
is closely related to another well-known problem. The cores of galaxy clusters are dense, 
and hence the gas there has very short cooling times compared to the Hubble time. This 
originally led to the idea of a cooling flow, and the expectation of large amounts of cold gas 
being deposited in such cores. However, recent Chandra and XMM-Newton spectroscopic 
observations (Tamura et al. 2001; Peterson et al. 2003) have found no evidence of X-ray gas 
cooling below ~ Tvir/3, suggesting a temperature floor of about 1 keV. Numerous scenarios 
have been envisaged in attempts to explain these observations. It has been argued that the 
gas may be cooling, but th e line emission is eith er absorbed, or greatly diminished through 
intrinsic X-ray absorption ( jPeterson et al.l 120011). or due to the presence of multiphase gas 
and inhomogeneity of iron abundance (IFabian et al.l l200ll ). An alternative to hiding the 
signatures of cooling gas is some form of non-gravitational heating. 

In order to offset the cooling in cluster cores, such heating must come from a cen- 
tral source. Burns (1990) found that 70 percent of cD galaxies contain radio sources at 
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as well as numerical simu- 



Active Galactic Nuclei (AGNs) can have a profound impact on their surroundings. Moreover, 
a num ber of tight correlations exist between the global galaxy properties s uch as spheroid 



mass (IMagorrian et al.l 1 19981 : iHaring fc Rixl |2004J ) and velocity dispersion (iGebhardt et al. 
2000), and those of the black holes found at their centres. 



AGN feedback provides a natural link between the black hole properties and those of the 
host galaxy, and has been invoked by a number of researchers to explain the disagreements 
between the hierarchical structure assembly paradigm and observations. Granato et al. 
(2004), Bower et al. (2006), Cattaneo et al. (2006) and Croton et al. (2006) have all 
investigated the role of feedback from radio galaxies in this process. These authors found 
that the feedback can suppress gas cooling and hence star formation significantly in massive 
hosts. However, the feedback models they used were largely phenomenological. Croton et al. 
assumed a rather arbitrary formulation for the rate of black hole fuelling, thus ensuring star 
formation is preferentially suppressed in massive galaxies. Cattaneo et al. similarly assumed 
that AGNs only switch on in haloes more massive than some critical value. Granato et al. 
assumed a constant fraction of AGN kinetic luminosity couples to the gas, and derived 
approximate kinetic luminosities from known scalings of observable quantities. Bower et al. 
postulated that cooling in galaxies is suppressed when the radiative output of the cooling 
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gas exceeds some fraction (typically ~ 0.5) of the Eddington luminosity. 



While encouraging, these approaches do not provide a realistic representation of the 
interaction between the radio jets and the surrounding gas. The aforementioned authors 
have to a large degree built into their models factors that guarantee the success of the AGN 
feedback process. Thus, whether or not radio galaxies really can provide sufficient heating to 
quench gas cooling by the required amount remains an open question. We address this issue 
by developing and applying a galaxy evolution model that, for the first time, incorporates 
detailed and physically motivated, radio source feedback. 

AGN feedback is an inherently intermittent process. Heating of the intracluster gas 
(ICM) decreases the rate of accretion onto the central black hole, until it is shut off completely 
(e.g. Shabala et al. 2008). Once the gas has had sufficient time to cool, the accretion can 
restart. It has been suggested (Jester 2005 and references therein) that accretion processes in 
Galactic X-ray black hole binaries (BHXBs) are analogous to those in AGNs, albeit on much 
shorter timescales. At least two different accretion states are observed in BHXBs; these 
are characterised by different spectral shapes and luminosities. The high luminosity /soft 
X-ray (thermal) spectrum state is identified with the standard Shakura & Sunyaev (1973) 
thin accretion disk. The resultant accretion flow is optically thick and radiatively efficient. 
By constrast, the low luminosity /hard X-ray spectrum state is identified with an optically 
thin accretion flow. The co oling times are long a nd trapped radiation is advected inward, 
resulting in a "puffed" disk (jPark fc Ostrikerll200ll ). The resulting quasi-spherical Advection 
Dominated Accretion Flow (ADAF; Narayan & Yi 1995) is thus radiatively inefficient, and 
the BHXB is observed in a low luminosity state. Radio jets are only observed in the ADAF 
phase, with quenc hing occurring when the BHXB enters the high/soft (i.e. thin disk) state 
jGallo et al.ll2003h . 



In this paper, we develop a galaxy evolution model that includes intermittent radio 
source feedback. In our model, existence of radio jets depends on the nature of the accretion 
flow solution, which in turn is sensitive to the interplay between cooling and heating of the 
cluster gas. We describe the major features of our galaxy formation model in Section [21 and 
incorporate feedback in Section [31 Constraints on the model are placed in Section [H and 
it is tested by means of a comparison with the statistical properties of the observed galaxy 
population in Section [51 

Throughout the paper, we adopt a flat cosmology of Q m = 0.25, fl\ = 0.75, h = 0.7 and 



as = 0.9, consistent with the 2dFGRS (ICoUess et al.ll200lh a nd WMAP fISeliak et al.ll2005h 
results, and the Millennium Simulation (ISpringel et al.ll2005l ). 
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2. Basic framework 

2.1. Dark matter halo evolution 

To follow the evolution of an individual dark matter halo, the analytical description 
of van den Bosch (2002) is employed. Van den Bosch uses an extended Press-Schechter 
formalism to derive halo mass accretion histories which are in good agreement with the 
ACDM simulations of Kauffmann et al. (1999a). The average mass accretion history of a 
halo with mass M in the local Universe is given by 



log 



M{z) 
M 



-0.301 



log(l + z) 

_log(l + Zi) 



where the fitting parameters z/f and z^ are a function of halo mass and cosmology. 
For the adopted flat cosmology, Zi is obtained by solving 



0.0055 



D{z,) 

Here, n{z) - "^^'+"^' 



fi^f + 0.283 [2((t2(0.254M) - ^^(M))] 



0.5 



(2) 



nA+nM(i+z)^ i S'^d D{z) is the linear growth factor normalised to unity at 
z = Q. Using the fit of Eisenstein & Hu (1999), 
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Following van den Bosch (2002), another fitting function is adopted for o"(M), 

/ (3.084 X lO-^F 
a(M) = as— ^ 







M 









/(32F) 

where F = 0.15 is the power spectrum shape parameter. The function / is 
f{u) = 64.087 X (1 + 1.074m°-^ - 1.581m°-^ + 0.954m°-^ - 0.185u° 



6N-10 



(3) 



(4) 



(5) 



Equations [2115] allow Zf to be determined, given a halo mass at redshift zero and cosmol- 
ogy. The second parameter describing the average halo mass accretion history is then given 
by 
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z/f = 1.211 + 1.858 log(l + Zf) + 0.308r]^ - 0.032 log 



M 



10"/i-iM 







(6) 



Together, equations [T}E] uniquely determine the average mass accretion history of a halo 
with a given mass at redshift zero. 

Detailed N-body simulations for the ACDM cosmology (e.g. Jenkins et al. 1998; Kauff- 
mann et al. 1999a) provide th e present-day halo rn ass function. Analytical fits to the results 
of the Millennium Simulation (jSpringel et al.ll2005l ) provide a halo mass function at each step 
of the evolution. Since our model does not explicitly follow mergers, such halo mass functions 
are required in order to predict galaxy population properties as a function of cosmic epoch. 



2.2. Galaxy formation and evolution 

2.2.1. Density profiles 

Gas will be accreted onto the halo along with the dark matter (e.g. White & Frenk 1991; 
Scannapieco et al. 2005). We take the mass of baryonic matter accreted at each point in the 
halo's evolution to be a fraction /b of the total accreted mass. For high mass haloes /b ~ 
Photoionization heating from the UV background reduces this fraction in low-mass haloes; 
this is discussed in Section 13.11 The gas is shocked to its virial temperature and density as 
it is accreted at the halo virial radius. 



Pvir 



T ■ 



3M„^ 



1 GM^i, fim}i 

2 Rvir 



1/3 



(7) 



where A,i, = ISn^ + 82{Q{z) - 1) - 39{Q{z) - \ f flBrvan fc Norman! Il998r ) is the mean 
overdensity at the virial radius; and /xmH is the mean particle mass, with \x = 0.62 for a fully 
ionised gas. The halo virial velocity is 



vir 



1/2 
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The potential well of a spherical dark matter halo is assumed to follow the Navarro, 
Frenk and White (NFW; Navarro, Frenk & White 1997) profile, 

Puuir) = - — (9) 



where scale factor Rg = for concentration index C = (^TIxW^M^^ ( Bullock et al. 



200ll ). Central density Pdm,o is set by the constraint that halo mass M^aio is contained within 



the virial radius. The NFW profile is amenable to integration in closed form, yielding 

Pdm,o = — , wi\ ' (10) 




where f{x) = [in (l + i) - j^] . 

We assume the density distribution of hot X-ray gas also follows Equation [TDl As with 
the dark matter, the constraint that all the gas is located within the virial radius determines 
the central density Pgas,o- In this fashion, the density profile can be determined at each point 
in the halo's evolution. 



2.2.2. Gas cooling and star formation 

The gas subsequently cools, with the cooling time of each gas parcel being given by 
= I fj,m^^^A(T z) ' ■^hsre electron density rig = ^^^^ for a fully ionised gas. This yields 

75/i mukBT ^^^^ 



2(2 + /i)2p...A(T,Z)- 

The cooling function A(T, Z) depends on gas temperature and metallicity. Detailed cooling 
models of Sutherland & Dopita (1993) for Z = O.IZq are used here. 

The cool gas is deposited onto an accretion disk after a dynamical time tdyn('") = y—, 
and subsequently star formation can take place. Kennicutt (1989) showed that star forma- 
tion rates in a sample of nearby spiral galaxies have a sharp cutoff below a certain surface 
density, in line with gravitational instability considerations. This critical surface density is 
representative of the disk gas density (except at the very edges), and thus one can define a 
critical disk mass below which no star formation takes place. Kauffmann (1996) gives the 
threshold density as S^nt = 0.59 ( ^^"^-i ) (^p^) ^opc"^. Integrating over the surface of 
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the disk, the critical mass is M„it = 7.5 x 10^ f^^^^^j (fj^) Adopting disk radius 
-Rdisk = 0.1 R^i^ (lKauffmannlll996l : ICroton et aLll2006l ). this yields 



= 7.5 X 10^ ' ^ 



0- 



^kms V V^pc 

Following Croton et al. (2006), the instantaneous star formation rate is taken to be 



(12) 



= asF- 



^dyn,disk 



(13) 



where disk dynamical time is tdyndiak = and star formation efficiency osf ~ 0.1 (Croton 
et al. 2006). Here, stellar mass refers to the total stellar mass within the dark matter halo. 



2.3. Black hole growth 

In standard quasar models the build up of black hole mass occurs via both accretion and 
mergers. In the present model the average accretion history of a given halo is considered, 
and thus mergers are implicitly included. When two progenitors merge, we assume the black 
holes also merge to form a new object with mass equal to the sum of the progenitor masses 
(e.g. Kauffmann & Haehnelt 2000; Croton et al. 2006). Thus the contribution from mergers 
to black hole mass build up follows the growth of the dark matter halo. 

Black holes will also grow through accretion of cold disk gas. Clearly, the growth rate 
will strongly depend on the availability of the cold gas. The shallower gravitational potential 
wells of lower-mass haloes make gas there more susceptible to ejection by supernovae and 
AGNs. Following Kauffmann & Haehnelt (2000; see also Malbon et al. 2007 and Lagos et al. 
2008) and parametrising the accreted gas fraction by halo mass, the mass accreted onto the 
black hole in a time At is 

l + (280kms-VK.)2' ^ ' 

3. Feedback processes 
3.1. Reionization 

At low halo masses, gas accretion is suppressed due to the ionizing radiation from 
extragalactic UV background photons. Gnedin (2000) modelled this effect by introducing 
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a filtering mass Mf below which the accreted gas fraction /b is reduced from the universal 
value /b,o = g^, 

Mz, M.,) = ^-^ —. (15) 

1 + 0.26 

Semi-analytic models typically assume accretion is suppressed when some measure of gas 
temperature (e.g. the average temperature, temperature of dense clumps or at the periphery 
of the halo) exceeds the virial temperature of the halo. Okamoto et al. (2008a) showed 
that such an approach overstates the effects of photoionization feedback, and it is necessary 
to follow the full merger history of the haloes. Since our model does not explicitly follow 
mergers, we adopt the Okamoto et al. (2008a) results for the filtering mass as a function of 
redshift. As is discussed in Section |U the effects of this photoionization heating are negligible 
at halo masses above 10^^ Mq. 



3.2. Supernovae feedback 

Newly formed massive stars end their lives as supernovae. The supernova outbursts in- 
ject a significant amount of energy into the surrounding gas, reheating the cold disk material 
and, in the case of the most powerful events, uplifting the gas out of the cluster. They also 
inject heavy elements into the gas. We follow the treatment of De Lucia et al. (2004) and 
Croton et al. (2006) in modelling these processes. 

At each point in time, a fraction of the stellar mass formed is instantaneously recycled 
into cold disk gas. The exact fraction depends on the assumed Initial Mass Function (IMF). 
Consistent with the adopted diet Salpeter IMF (see Section 14.11) we take /^^c = 0.33. 

The amount of reheated disk gas is proportional to the recycled mass, and hence to the 
mass of stars formed in the timestep of interest, 

AM,,._,,d = Cdi^kAM,, (16) 



where e^isk = 3.5 (ICroton et al.ll2006l ). This raises the total thermal energy of the hot gas by 

AE,„, = ^AM,„,Kl- (17) 
The energy output of a supernova event is 

A^sN = ^e^.,AM,V,l, (18) 
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where T^n = 960 kms~"^ is found by equating |V^^ with the mean energy in supernovae 
ejecta per unit mass formed, assuming a diet Salpeter IMF and an energy output of 10^^ J 
per supernova event. The efficiency ehaio with which disk gas is reheated is a free parameter 
that is constrained by observations. If this energy exceeds that required for the heating, a 
fraction of cold gas will be ejected outside the halo, 

vir 

The ejected gas is returned within a few halo dynamical times, 

A Ti -^'^cjcctcd^^ /'0^^ 

^-^•^ returned ^returned T ; V"^^/ 

dyn,halo 

where tdv„haio — and e.^tumcd ~ 0.5. If the available energy is less than that required for 
heating, no gas is ejected, and the reheated mass is 

AM„ = (21) 

vir 

It is worth noting that the exact values for most parameters used in this implementation are 
not important, as the supernovae feedback mechanism is parametrised by the quantity ehaio, 
the magnitude of which is constrained by observations. 



3.3. Radio source feedback 



Apart from forming stars, a fraction of the cold gas will be accreted onto a disk in 
the vicinity of the black hole. In order for the accreted gas to spiral in towards the centre, 
conservation of momentum requires its angular momentum be removed. Although magnetic 
fields could facihtate this process, in practice this requires ejection of gas. The power of the 
resultant jets depends on the amount of fuel available. When powerful enough, these jets 
will terminate in hotspots, and backflow of radio plasma will inflate the cocoons observed in 
powerful FR-II radio sources. The radio source will then expand, driving a shock through 
the intraclustcr gas and thereby heating the swcpt-up gas located between the contact dis- 
continuity delineating the cocoon and the bow-shock. This heating will in turn limit the rate 
of gas cooling and subsequently the rate at which the black hole is fuelled, until production 
of jets of sufficient power to affect the gas is shut off. Once the gas has time to cool, the 
process is restarted. Below we develop a detailed model for following the interplay between 
the heating and cooling of the gas. 
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3.3.1. Shock heating 



Observations (e.g. Leahy et al. 1989; Subrahmanyan et al. 1996) suggest radio sources 
are self-similar, meaning the bow-shock and cocoon radii are related by -Rcocoon = ^Rshock- 
We adopt the models of Kaiser & Alexander (1997) and Alexander (2000, 2006) to describe 
the evolution of the cocoon radius -Rcocoon with time. The gas density profile in each halo 
is approximated as a flat core — double power law atmosphere. In other words, the source 
initially evolves in a flat atmosphere, followed by two power law profiles of the form p(r) = 

. This would correspond to an expansion through a high-density core, into the 



galaxy, and then the cluster gas. 



Both analytical work (e.g. Kaiser & Alexander 1997) and numerical simulations ([Reynolds et al 



200ll ) suggest the swept-up gas lying between the shock and cocoon radii is isobaric. Such 
an arrangement is facilitated by backfiow of swept-up gas from the hotspot to the sides of 
the cocoon. We therefore follow Heinz et al. (1998) and Alexander (2002) and assume the 
cocoon to be spherical. 

We further assume the swept-up gas evolves isothermally, consistent with cluster obser- 
vations. Following the treatment of Alexander (2002), the supersonic expansion of the radio 
source shocks this gas to a temperature 



15 3 - 



3-/3 / ^imn \ ]_^2 



16 11 - /3 V y A2 



cocoon ' 



(22) 



15 



where A and /3 are related by 1 — A'^ — 4(ii„^) • 

Important parameters governing the radio source expansion and hence the heating of 
the gas are the gas density profile, which is given by Equation [HI cocoon axial ratio R^, which 
is r elated to the jet opening angle and is fixed at -Rt = 2, a value appropriate for Cygnus 
A (IKaiser fc Alexanderl 119971 : iBegelman fc Cioffil Il989l ): and jet power Qj^t- Following the 
dynamical model of Kaiser & Alexander (1997), the cocoon radius is 



■ ^ N 3/(5-/3) 

-^cocoon ('^) ^D-^corc ( 1 



(23) 



where r = f ^^:£2^£^£is.\ jg ^ convenient timescale; and the dimensionless constant a^i is 

\ ''je* 1 II 1 

given by (IKaiser fc Alexanden 119971 ) 



an 



(r. + i)(re-i)(5-/9)3 

187r (9 [F, + (F, - l)i?2 ] _ 4 _ 



1/(5-/3) 



(24) 
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The adiabatic indices for the cocoon and external gas, Tc and Fx respectively, are equal to 
5/3, corresponding to non-relativistic material. 



3.3.2. Jet power 

Theoretical work (Narayan et al. 1998; Meier 2001) as well as observations of black hole 
X-ray binaries (Kording et al. 2006; Fender et al. 2004) suggest that in both these objects 
and AGNs at least two different accretion states exist. In both cases angular momentum 
is removed from the accretion disk by viscosity, allowing the gas to spiral in towards the 
central black hole. At high inflow rates (com pared to the Eddington rat e), the accretion flow 
is described by a standard thin disk solution (IShakura fc Sunyaevlll973l ). Here, the accretion 
disk is geometrically thin and optically thick, and produces a quasi-blackbody spectrum. 
This flow is rad iatively efficien t, as all the energy released through viscous dissipation can be 
radiated away (INarayanl |2002| ) . By contrast, at low accretion rates the flow is geometrically 
thick and optically thin. As a result, the cooling times are long, and instead of being radiated 
away the thermal energy of the inflowing gas is advected inward. The resultant Advection 
Dominated Accretion Flow (ADAF; Narayan & Yi 1995) is thus radiatively inefficient, and 
can produce powerful jets. 

Meier (2001) gives the jet power generated by a geometrically thick, optically thin ADAF 
disk as 



Q 



jot, ADAF 



1.3 X 10 



^adaf(0-55/adaf + 1-5/adafJbh + Jbf 



W 




(25) 



for rhj: 



< m^rit, where merit is a parameter discussed in the following section. 



Here, Mbh is the black hole accretion rate, and M„ 



2.3 

^rad 



Meyr- 



fradC^ frad V IO^Mq , 

is the Eddington accretion rate for a radiative efficiency e^ad- Disk viscosity is responsible 
for transporting angular momentuna outward, and is given by the coefficient Oadaf ~ 0.3 
( iNarayan et al.l Il998l : INarayanl |2002| ). Coefficients (7adaf and /adaf refer to the ratios of 
actual angular velocity and azimuthal magnetic field to those calculated by Narayan & Yi 
(1995), and Jbh is the black hole spin. Following Meier (2001) and Okamoto et al. (2008b) 
we adopt (^adaf = 2.3, /adaf = 1 and Jbh = 0.5, noting that the jet power is not very 
sensitive to realistic changes in these values. For this choice of parameters. Equation [25] can 
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be rewritten as 



Q. 



jct.ADAF 



4.5 X 10^9 ' 



Moyr 



-1 



W 



0.79Mr„c2 



(26) 



For accretion rates Mbh > "/Ticrit^odd the accretion flow is described by the standard Shakura- 
Sunyaev thin disk solution, which for the adopted viscosity parameter is given by 

0.2 



Q = 6.0 X 10"^ ( 



-0.3 



Moyr-i 



(27) 



This mode of accretion (referred to as the "quasar mode" by Croton et al. 2006) is radiatively 
efficient, and typically result s in weaker lets than the ADAF cas e. Thus it is identified 
with a luminous AGN phase (jCroton et al.ll2006l : iBower et al.ll2006h : while the ADAF phase 
corresponds to the "radio mode" of Croton et al. and is responsible for radio source feedback. 
In what follows, we assume that when the instantaneous black hole accretion rate exceeds 
the critical rate, the disk is in the "quasar mode"; and when the instantaneous accretion 
rate is below this value, a radio jet with power given by Equation [26] is generated. We 
further assume that in this case the actual accretion rate onto the black hole is given by 
Mbh = ""T-crit^odd, i-e. the disk is fuelled at the maximum rate (e.g. Hopkins et al. 2006). 
The usual radiative efficiency of e^^d = 0.1 is assumed in evaluating M^dd- 



3.3.3. Radio source intermittency 



Croton et al. (2006) and Bower et al. (2006) both took the "radio mode" heating 
provided by the jets to be continuous. Observations of black hole binaries and the apparent 
balance between heating and cooling rates, however, suggest otherwise. Jets will be inter- 
mittent for two reasons. Feedback limits the availability of cold gas required to fuel the 



central engine ( jShabala et al.ll2008l ). resulting in weaker jets that are more liable to disrup- 
tion by instabilities in the dense central regions. On the other hand, when the radio cocoon 
expands outside the high density regions, the central gas cools rapidly. The resulting high 
black hole accretion rates qualitatively change the nature of the accretion flow from a ra- 



diativ ely inefficient ADAF to a radiatively efficient standard thin disk flow (iNarayan et al. 
19981 ). This once again results in a much weaker jet (by ~ 3 orders of magnitude according 
to Equations [2S] and [27j) which can be disrupted. 



In the assumed density proflle, the jets are most likely to only be susceptible to disrup- 
tion by Kelvin-Helmholtz instabilities in the denser central regions of the halo (lAlexander 



-14- 



20021 : iKaiser &: Alexanderl 119971 ). Following the analysis of Alexander (2002), the minimum 
power a jet must have to traverse the fiat core is taken as 



Q 



jct,min 



10 



37 



Rr 



/-^corc 



2.5 kpc/ \ 1.5 X 10^2 i^gjn- 



W. 



(2^ 



Similarly, for a radio source that has expanded outside -Rcoro, the jet power must exceed 
Qjct,min oc -Reo~fonPcorc in ordcr to reach the hotspots. If this does not happen, or the instan- 
taneous rate of accretion onto the black hole men exceeds rh^^it causing the disk to become 
radiatively inefficient, the jet is terminated, and a new jet inflates a new radio cocoon. For 
a radio source that is still within the core, i.e. when -Rcocoon < -Rcoro, the minimum required 
jet power is given by replacing the core radius with the cocoon radius in Equation [281 

The flnal important point concerning jet intermittency relates to the critical accretion 
rate denoting the transition between radiatively efficient and inefficient disks. Narayan et al. 
(1998) estimate that an ADAF will switch to a radiatively efficient (thin disk) flow when the 
accretion rate in Eddington units exceeds m, 



critjUp 



a 



ADAF 



0.1 for the adopted viscosity 
value. This threshold is consistent with the observed transition from low/hard to high/soft 
state in X-ray binaries (lEsin et al.l 1 19981 : iKording et al.l 120061). Intere stingly, the reverse 
transition is less well deflned, with rhcnt.down ~ 0.01 — 0.1 ( lNarayarul2002l ). The thin disk will 
develop a hole in its inner regions, where an ADAF is allowed. As the accretion rate drops 
further, the region hosting an ADAF will grow, until this becomes the dominant mode of 
accretion. Hence, in the present model we treat d^^^ as a free parameter in the range 
0.01-0.1. 



4. Constraints on model parameters 

There are three free parameters in the presented model. Supernovae heating is parametrised 
by e^aio, the efficiency with which supernovae heat the cold disk gas. The other two parame- 
ters relate to AGN feedback. Fraction of cold disk gas accreted onto the black hole is given 
by e^,,,;., and m„itd„„„ denotes the dimensionless accretion rate (in Eddington units) below 
which the accretion flow switches from a thin disk to an ADAF solution. 

Constraints on these parameters come from observations. Radio source heating is only 
important in massive haloes and at late times, where the gas cooling rate can provide suf- 
flcient fuel for a radio source. This can be seen in Figure [T], where evolution of the stellar 
component is shown for haloes with masses 10^^ and 10^^ Mq at the present epoch. In the 
10^^ Mq halo, supernovae feedback and reionization heating result in the stellar content of 
the halo being at all times consistently lower by a similar factor than when only gas cooling 
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is included. By contrast, this has httle effect on the 10^^ Mq halo, and even then only at the 
highest [z > 3) redshifts. However, at z ~ 1 the cooling rate becomes high enough for radio 
source feedback to become important. From that point it is AGN feedback that dominates 
the cooling history of the halo. 




12345 012345 



(a) 1.6 X 10" M0 (b) 1.6 X lO^^ Mq 

Fig. 1. — Evolution of the stellar mass content in haloes with redshift zero masses of 1.6 x 10^^ 
and 1.6 X 10^^ Mq. Different curves correspond to different models: no feedback (short- 
dashed); reionization and supernovae heating (long-dashed); and full model including AGN 
feedback (solid line). Halo mass build-up (dotted line) is plotted for reference. All curves are 
for e^aio = 0.01, e^cc = 0.02, mcrit.down = 0.03. Reionization and SNe feedback are important 
at early times and in low-mass haloes, while radio source feedback dominates at late epochs 
in massive hosts. Thus, on scales of interest the reionization and supernovae heating (long- 
dashed) curve coincides with the full model (solid) curve for the 1.6 x 10^^ Mq halo; and the 
no feedback (short-dashed) and reionization plus supernovae (long-dashed) curves coincide 
for the 1.6 x 10^^ M© halo. 

Contributions of individual feedback components to halting the build-up of stellar mass 
are shown in Figure [H Feedback from supernovae and reionization is dominant in low-mass 
(< 10^^ Mq) haloes, while star formation is significantly impaired by AGN feedback in 
massive haloes. 



4.1. Initial Mass Function 



Both the rest-frame UV flux used to derive star-formation rates, and opti cal and infra- 
red emission tracking the integrated stellar mas s are heavily obscured by dust (ISteidel et al. 



19991 : iMadau et al.lll998l : iGiavalisco et al.ll2004j ). Fortunately, the rest-frame near infra-red 



iT-band luminosity function is almost a magnitude less sensitive to dust extinction than at 
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Mhalo/Msun 



Fig. 2. — Stellar mass at the present epoch for various models: cooling (dotted); reionization 
(short-dashed); reionization and supernovae (long-dashed); full feedback model including 
AGN (sohd line). The curves are again for e^aio = 0.01, e^^^ = 0.02, m^^i.tiown = 0.03. AGN 
feedback is important in massive (M^aio > 10^^ M©) galaxies. 



other wavelengths (IKauffmann et al.lll999al ). making it a very useful probe of stellar mass. 



The main uncertainty in converting between stellar masses and luminosities therefore 
comes from the adopted IMF. The standard relation frequently used is the Salpeter (1955) 
IMF given by ^^^^^^ oc M~^'^^ derived for solar neighbourhood stars. It has recently been 
argued, h owever, tha t field stars do not provide a good proxy for the low-mass end of 
the IMF (jScalolll998l ) due to the sensitivity of the derived mass function to star formation 
histories at these m asses, as well as being at odds with H a luminosities of galaxy populations 
(IBaldry et al.ll2004l ). Studies of star cluste rs and gala xies suggest that a "universal" IMF 
(in a statistical sense) may be applicable (jScald Il998l ). To this end, a number of models 
have been put forward. Kennicutt (1983), Kroupa (2001) and Chabrier (2003) IMFs are in 
good agreement with the Salpeter one at high masses, but predict less low-mass stars. Bell 
et al. (2003) adopt a "diet Salpeter" IMF, which has the same colours and luminosities as 
the normal Salpeter IMF, but only 70% of the mass due to a lower number of low-mass 
(and therefore faint) stars in order to match observed rotation velocities in the Ursa Major 
cluster. Using the vertical velocity dispersion of stars in disk-dominated galaxies, Bottema 
(1997) argued for an ev en lower mass-to-lig ht ratio. As a result, the uncertainty in the IMF 
slope is as much as 0.5 (IBaldry et al.l 12004 ). 
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Bruzual & Chariot (2003) compared the spectral properties of a number of IMFs, and 
found that the [B — V) and ( V — K) colours predicted by these are very similar. The main 
differences arise in the mass-to-light ratios. Bell et al. (2003) estimate that the choice of 
IMF results in the mass-to-light ratio varying by up to 0.5 dex, with top-heavy IMFs having 
lower mass-to-light ratios. 



4.2. Efficiency of supernovae feedback 

4-2.1. Counts of faint galaxies 

As Figures [1] and [2] show, low-mass haloes are little affected by AGN feedback. Since 
these haloes host galaxies with relatively low stellar masses, the low-mass end of the stellar 
mass function (or, correspondingly, the faint end of the iC-band luminosity function) can be 
used to constrain the contribution of supernovae to suppressing star formation. 

In Figure[3]the local stellar mass function is plotted for various values of Chaio- The curves 
are constructed by convolving the star formation histories of individual haloes with the halo 
mass function obtained at each redshift from the Millennium Simulation (Section l2.1l) . Shown 
by open symbols are the stellar mass functions derived from the Cole et al. (2001) 2dFGRS 
iT-band luminosity function assuming Salpeter and Bottema IMFs. These correspond to the 
upper and lower limits, respect ively, on the me an mass-to-light ratio used to convert from 



luminosities to stellar masses (IBell et al.l l2003l ). In the case of Salpeter IMF, we use the 



average mass-to-light ratio of 1.32 as quoted by Cole et al. For the Bottema IMF this value 
is reduced to 0.82. The filled symbols correspond to the stellar mass function derived with 
a die t Salpeter IMF, obtained by decreasing the Salpeter mass-to-light ratio by 0.15 dex to 



1.17 (IBell et al.ll2003l ): this gives good agreement with the Cole et al. stellar mass function 



derived from J-band observations. 

It can be seen that e^aio = 0.01 provides a good match to the local stellar mass function 
at the low-mass end. Although the choice of IMF impacts the derived mass function at 
the bright end, in the absence of AGN feedback the number densities of the most massive 
galaxies are overpredicted significantly. 



4-2.2. Distribution of stellar mass 

Low mass haloes are much more abundant than their higher mass counterparts- for 
example, 10^^ M© haloes are five orders of magnitude more abundant than 10^^ Mq ones 



J I I I I I I I I I I I I I I I I I I I I I I l_ 

'9.5 10 10.5 11 11.5 12 
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Fig. 3. — Local stellar mass function for various supernovae feedback efficiencies. Observed 



mass functions are derived from the 2dFGRS J-band near infra-red observations ( Cole et al. 



20011 ). assuming a Salpeter (open circles), Bottema (open triangles) or diet Salpeter (filled 
squares) IMF. Model curves are: cooling (dotted), reionization (short-dashed), reionization 
plus supernovae feedback with ehaio = 0.01 (solid), e^aio = 0.04 (long-dashed). Bright end 
counts are significantly overpredicted without AGN feedback. A colour version of this Figure 
is available in the online article. 
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( iJenkins et al.ll200ll ). Since massive haloes host the more massive galaxies, however, the 
question of how the total stellar mass is distributed between haloes is a non-trivial one. We 
plot in Figure H] the stellar mass fraction as a function of halo mass at the present epoch. 
Haloes with masses less than 3 x 10^^ — 3 x 10^^ Mq are seen to dominate the total stellar 
mass budget, despite the total mass distribution peaking at higher halo masses. 
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Fig. 4. — Fraction of stellar mass contained in different mass haloes at the present epoch. 
Dashed curve represents the cooling-only model. Solid curve includes reionization heating, 
supernovae feedback for e^aio = 0.01, and AGN feedback with e^^c = 0.02 and = 0.03. 

Fractional distribution of halo mass (dotted curve) is plotted for comparison. 



This result suggests that the star formation history of the Universe is largely determined 
by reionization and supernovae, rather than AGN, feedback. This is certainly true for z > 1, 
when AGN feedback does not greatly affect the evolution of the stellar mass content even 
in the most massive haloes (Figured]). Hence a guide to the efficiency of the SNe feedback 
can be obtained from observations of the evolution of the star formation rate (SFR) density, 
plotted in Figure [51 Observations at 2; > 1 suggest e^aio ~ 0.01. It is also clear from Figure [5] 
that models without AGN feedback significantly overpredict star formation rates at later 
epochs. The observed star formation rates are derived from a combination of UV, IR, O [II], 
radio and X-ray observations (see caption) assuming a diet Salpeter IMF. It is worth noting 
that, due to the assumed IMF, these are close to the upper limits on the actual star-formation 
rates at each epoch, and the actual SFRs can be up to 0.35 dex lower than these values. 
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Fig. 5. — Evolution of the global star formation rate density. Various curves denote models 
with: cooling only (dotted); reionization (short-dashed); and reionization plus supernovae 
with e^aio = 0.01 (solid) and e^aio = 0.04 (long-dashed). Observations are from the compilation 
of Hopkins (2004) and J-band observations of Cole et al. (2001; open squares), assuming a 
diet Salpeter IMF. Hopkins (2004) points represent various UV (open circles), radio (filled 
squares) and X-ray (filled circles) observations. A colour version of this Figure is available 
in the online article. 



4.3. AGN fuelling efficiency 



With supernovae feedback efficiency fixed, the fraction of gas accreted onto the central 
black hole can be determined by invoking known observed relations between the mass of the 
central black hole and global halo and galaxy properties. One of these is the correlation 
between the local masses of the black hole and stars in the bulge (as derived from their 
luminosities; Magorrian et al. 1998; Haring & Rix 2004). 




Fig. 6. — Mbh — Mbuige relation in the local Universe as a function of accretion efficiency. 
Curves are: e^^^ = 0.01 (short-dashed); e^^^ = 0.02 (solid); and e^^^ = 0.05 (long-dashed line). 
All models include reionization and supernovae feedback with e^aio = 0.01, but no AGN 
feedback. Data points are from Haring & Rix (2004). 



Although an estimat e of the stellar content in the bulge can be obtained by invoking 
disk stability arguments (lEfstathiou et al.lll982l ). such a formulation ignores mergers and 
thus attributes too much stellar mass to the disk. Benson et al. (2007) determined spheroid- 
to-disk light ratios in a sample of ~ 9000 SDSS galaxies, and found that ~ 0.3 - 1.0 
for galaxies with stellar masses between 10^ and 10^^ Mq. We adopt the Benson et al. 
relation between stellar mass and the fraction of stars located in the spheroid. Figure [H] plots 
the Mbh — Mbuigc relation thus determined for models without any AGN feedback. Accreted 



fractions 
galaxies. 



0.02 — 0.05 are required to match observations of low-mass (Mtuig,, < 10^^ 



Because our model does not account for galaxy groups and clusters, here the model 
bulge masses refer to the sum of stellar masses of all the galaxies in a given halo, multiplied 



-22- 



by the bulge fraction. Although this fraction is likely to vary between objects (in particular, 
satellites are more likely to be disk-dominated than cD galaxies), one would expect such 
bulge estimates to be reasonably accurate because of the dominant contribution of the central 
galaxy to the stellar mass budget. 

Due to the self-regulatory nature of the feedback processes, even in the presence of AGN 
feedback the nature of the M^^ — M^uige relation at the high-mass end is not expected to 
be altered substantially. Rather, one would expect both the black hole and bulge mass to 
be reduced by a similar factor relative to the no-feedback case. This point is considered in 
Section 15. 1[ 



4.4. Thin disk — ADAF transition 

The dimensionless accretion rate m„it regulates how often jets are re-triggered. An 
observational constraint on the efficiency of feedback is provided by the bright end of the 
optical luminosity (or stellar mass) function. 
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Fig. 7. — Local stellar mass function for various mcrit,down values. Reionization and supernovae 
feedback for e^aio = 0.01 are included. Curves are for: m„it,down = 0.03 (solid), mcrit,down = 0.07 
(dashed). The no- AGN prediction (dotted curve) is plotted for comparison. Observed mass 
functions are as in Figure [31 A colour version of this Figure is available in the online article. 

Figure [7] shows the z = stellar mass function for a diet Salpeter IMF for a range of 
'^crit,down values. The best fit is given by m„it,down ~ 0.4. Importantly, this set of plots shows 
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that realistic AGN feedback can resolve the number counts problem at the high-mass end of 
the stellar mass function by suppressing star formation in massive haloes. 

4.5. AGN feedback models 

A number of authors (e.g. Croton et al. 2006; Bower et al. 2006; Cattaneo et al. 2006; 
Granato et al. 2004) have implemented AGN feedback in galaxy formation models. The 
major difference between their results and ours lies in the intermittency of our physically- 
motivated prescription. Figure M compares our predictions to the model of Croton et al. 
(2006). It is worth noting that the value of the scaling constant describing the strength of 
AGN feedback (kagn in Croton et al.'s Equation 10) is different to the value obtained by 
those authors. This is a result of the differing model implementations. 




Fig. 8. — Comparison of AGN feedback models. Left panel: Stellar masses at 2; = 0. Solid 
curve represents our model with the same parameters as in Figure [2l Results from Wang 
et al. (2006) from the Millennium Simulation, with feedback as implemented by Croton et al. 
(2006), are shown by a dotted curve. Middle panel: Cooling histories for a 10^^ halo. 
Right panel: Accretion rate onto the black hole, in Eddington units. We predict lower rates 
of accretion than the model with continuous feedback and decreases at lower redshifts. 
The epoch of powerful AGN activity begins at 2; ~ 1.4 as man falls below merit, down = 0.03. 

Since both models reproduce the observed local stellar mass function, the stellar masses 
for each halo at the present epoch are in agreement (left panel). Cooling histories (middle 
panel) are also similar. However, the two models differ in their predictions of the rate at 
which cold gas is accreted onto the central black hole (right panel). The Croton et al. model 
predicts an increasing (in Eddington units) accretion rate; while our model predicts men to 
decrease with time. At z ~ 1.4 the accretion rate falls for the first time below m„,it,down, 
triggering a powerful jet episode and resultant intermittent activity. The interplay between 
gas heating, cooling and ejection by powerful AGN outbursts is crucial to explaining both the 
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observed duty cycle of radio galaxies and quasars (Shabala et al. 2008) and their evolution. 
We defer a detailed discussion of this point to a future paper. 

5. Cosmic downsizing and suppression of star formation 

5.1. Mbh - Mb^ige relation 

Figure [7] shows that the model can explain the observed dearth of massive galaxies and 
match the stellar mass function at the present epoch. It must also reproduce the observed 
local correlation between black hole mass and global properties such as spheroid mass. This 
is plotted in Figure M 




Fig. 9. — Mbh— Mbuigj, relation at z=0 for the best-fit (eh^io = 0.01, e^^ = 0.02 and m„it = 0.03) 
model (solid line). Data points are again from Haring & Rix (2004). Prediction for the no- 
AGN case (dashed line) is plotted for comparison; this corresponds to the solid curve in 
Figure [61 

Although the final stellar and black hole masses in each halo are lower than when 
feedback is not included, these clearly grow "in step", with the slope of the final relation 
remaining unchanged. Thus radio source heating provides a self-regulatory mechanism that 
can simultaneously limit the growth of the central black hole and the bulge by a similar 
amount. 
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5.2. Evolution of global properties 

5.2.1. Star formation rate 

Figure dU] plots the star formation rate density as a function of cosmic epoch for the 
model providing the best fit to the local stellar mass function (Figure [7]), e^aio = 0.01, 
Cacc = 0.02, m„it down = 0.03. The no-AGN predictions are also plotted for comparison. As 
argued in Section 14.2.21 models without AGN feedback significantly overpredict the global 
star formation rate density at the present epoch. Inclusion of such feedback, on the other 
hand, provides results that are in good agreement with observations. 




0.5 1 1.5 2 2.5 3 3.5 4 4.5 
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Fig. 10. — Evolution of the global star formation rate density. Solid curve represents the 
best-fit model {e^-^^i^ = 0.01, e^cc = 0.02, m^rit.down = 0.03) as given by the local stellar mass 
function. Dotted curve is for the same parameters but without AGN feedback. Observed 
points are as in Figure O 



Integrating the predicted local stellar mass function yields a total present-day stellar 
mass density of = 4.5 Mq Mpc~^, or = 3.3 x 10~^. This is consist ent with observational 
estimates of = (2.8 ± 0.4) x 10^^ Jcole et al.lboOll : ISell et al.l[2003h and (3.5 ± 0.4) x lO'^ 
(iKochanek et al.ll200lh made with the same diet Salpeter IMF that is used in the models. The 
integral of the predicted star formation rate density over the Hubble time is fi^ = 5.1 x 10^^, 
since a fraction /.^^ = 0.33 of stars terminate in supernovae. 
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5.2.2. Stellar mass function 



A successful galaxy formation model must correctly describe the evolution of the galactic 
population, while simultaneously matching the local observable properties such as the mass 
function. This has traditionally proved very difficult, with most models either reproducing 
the large number of massive galaxies at high redshift but significantly overpredicting the 
bright end of the local luminosity function (e.g. Kauffmann et al. 1999a,b); or conversely 
matching the local K-band luminosity function but underpredicting the bright galaxy counts 
at z> 0.5 (e.g. Cole et al. 2000; Baugh et al. 2005). 

Figure [TT] plots the evolution of the stellar mass functi on predicted by the m odel. It is 
in remarkable agreement with the results of the K20 su rvev (IFontana et al.ll2004f) . and those 
of the FORS Deep Field (FDF) and GOODS South jDrorv et al.l l2005h survevs. out to a 
redshift of z > 1.5. The only other model in the literature to date providing a similar match 
is that of Bower et al. (2006). 

The last two panels of Figure [11] illustrate that the model consistently underpredicts 
galaxy counts at 2; > 2. The uncertainty in the adopted universal IMF cannot account for 
this as there is a good match between the model and observations at low redshift. However, 
observations of sub-mm and Lyman break galaxies (IBaugh et al.l |2005| ) as well as element 
abundances (INagashima et al.l 120051 ) suggest that a more top-heavy IMF is required at 2; > 
2, yielding lower mass-to-light ratios and hence rendering the derived high-redshift stellar 
masses as overestimates. The physical motivation for two different IMFs comes from a 
change in the dominant star formation mechanism with cosmic epoch. At high redshifts star 
formation is mainly trig gered by mergers; while at the pres ent epoch it follows the standard 
disk formation scenario (lKauffmannlll996l : iKennicuttI 1 1 9981 ) . 



Figure [12] shows the evolution of the stellar mass function under the assumption of 
a Bottema IMF. This gives good agreement at high redshift, but overpredicts the stellar 
content at low z. These results suggest a switch from the diet Salpeter IMF to a more top- 
heavy Bottema IMF at z > 2. Baugh et al. (2005) estimate that quiescent star formation 
(i.e. in disks) dominates for z < 2 — 3, consistent with our findings. 

The slight overprediction in the counts of low-mass galaxies at 2; > is consistent with 
the results of other models (e.g. Bower et al. 2006). There are two possible causes of 
this. One possibility would be current implementations of supernovae feedback truncating 
star formation too early. Alternatively, observational constraints at the low-mass end are 
likely to suffer from incompleteness, as evidenced by substantial field-to-field variations in 
Figures [11] and [12] 
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Fig. 11. — Evolution of the stellar mass function (solid lines). The corresponding quantity 
at 2; = is shown as a dotted line for comparison. Observed data points are derived from 
multiwavelength photometry and use of Bruzual & Chariot (2003) spectral synthesis models 
with a diet Salpeter IMF. Filled triangles are from Cole et al. (2001); filled squares and 
circles represent Drory et al. (2005) GOODS-South and FORS Deep Field observations; 
while the K20 survey of Fontana et al. (2004) is shown by open symbols. 




Fig. 12. — Same as Figure [TT] but for a Bottema IMF. A change to a top-heavy IMF is 
required for z > 2 in order to match observations. 
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5.3. Cosmic downsizing 



From Figures [TT] and [12] it is clear ttiat tlie liigli-mass end of tlie local stellar mass 
function is assembled at 2; > 1, while lower- mass objects continue forming stars until the 
present day. This is consis tent with observations of the most massive galaxies having red 
colours (jPanter et al.l 120071 ) and in stark contrast to the hierarchical assembly of halo mass. 



The left panel of Figure [13] shows the evolution of star formation rate density as a 
function of stellar mass. The stellar mass bins are chosen to follow the observations of 
Juneau et al. (2005) and Brinchmann et al. (2004) converted to a diet Salpeter IMF. The 
same quantities for a model with no AGN feedback are plotted in the right panel. 





(a) AGN feedback 



(b) No AGN 



Fig. 13. — Star formation rate densities broken up by stellar mass for the best-fit (left) 
and no-AGN (right panel) models. All curves are for: 10^'^^ < ^ < 10^°-^^ (solid line, 
open triangles), 10^°-^^ < ^ < 10"'°^ (dashed line, closed circles), 10"-°^ < ft < ^O"'^^ 
(dotted line, open squares). Data points for z > 0.5 are from Juneau et al. (2005); and from 
Brinchmann et al. (2004) for the local volume. A colour version of this Figure is available 
in the online article. 



Figure [12] shows that inclusion of AGN feedback changes the galaxy assembly paradigm 
from hierarchical to anti-hierarchical. There are three related reasons for the preferential 
suppression of star formation in massive galaxies at low redshift. Firstly, lower densities in 
the low-z Universe mean the ratio of cooling time to free-fall time is higher there than at high 
z, making AGN feedback more effective at low redshifts. Furthermore, black holes are more 
massive at low z, allowing higher jet powers to be generated by ADAF disks (Equation [26]) . 
Apart from providing greater shock heating, this also introduces an ad ditional mechanism of 
star formation suppression via gas ejection by powerful radio sources (IBasson fc Alexander 



2003]). The absence of this mode of feedback in the presented model is most likely responsible 
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for the disagreement between model predictions and observations in intermediate and massive 
galaxies and at low z. 
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Fig. 14. — Mean stellar ages weighted by mass (left) and luminosity (right panel) for models 
with reionization and supernovae feedback but no AGN (dashed); and the full model includ- 
ing radio source heating (solid curve). When AGN feedback is included the population is 
bimodal, with the most massive galaxies containing oldest stars. 



A final illustration of radio source heating giving rise to a massive, red population of 
galaxies is given in Figure [TH Here the mean mass- and luminosity-weighted stellar ages at 
the present epoch are plotted as a function of galaxy stellar mass. Results for a model with 
no AGN, but including gas cooling, and reionization and supernovae feedback are also shown 
for comparison. In the absence of any feedback the most massive galaxies host the youngest 
stars, and hence are blue, as expected from the standard hierarchical formation models. The 
situation changes dramatically when radio source heating is introduced, with galaxies more 
massive than 10^^ Mq at 2; = containing older stars than predicted by the no-AGN case. 
Such bimodality is consistent with the semianalytical model results of Croton et al. (2006) 
and Bower et al. (2006) and, significantly, the observed sha rp break between old an d young 
stellar populations at ~ 3 x 10^° Mq in SDSS galaxies (IKauffmann et al.ll2003bl ). 



6. Summary 

In this work, we present a new galaxy formation and evolution model. Gas is assumed to 
be accreted onto the halo along with the dark matter, and distributed in a quasi-isothermal 
fashion. Gas cooling competes with feedback heating from the ionizing UV background, 
supernovae and radio sources. Low-mass haloes dominate the stellar mass budget of the 
Universe at the present epoch, and efficiency of supernovae feedback is tightly constrained 
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by the star formation history of the Universe and counts of low-to-intermediate-mass galaxies. 
The observed local correlation between central black hole masses and those of stars in the 
bulge sets the rate at which gas is accreted by the black hole. 

For the first time, we implement realistic intermittent AGN feedback resulting from 
shock heating by expanding radio sources. Radio source heating is found to suppress star 
formation in the most massive galaxies, providing an excellent match to the observed local 
stellar mass function. Observations of the stellar mass function at earlier epochs suggest 
the stellar IMF may change qualitatively around z = 2 from a merger-dominated, top-heavy 
IMF to one consistent with observations of the local Universe. 

The model reproduces the anti-hierarchical paradigm of structure formation, with the 
top end of the stellar mass function in place by z ~ 1.5. Also predicted is the sharp observed 
transition in the mean stellar ages of galaxies at Mj, ~ 2 — 3 x lO"*^" Mq, with the more 
massive structures being "red and dead" . The self-regulatory nature of this AGN feedback 
channel is illustrated by our model reproducing the local black hole - spheroid mass relation, 
and suggests that intermittent radio source feedback can provide the missing link in galaxy 
formation and evolution models. 
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